Data relationships revisited

I tried transforming the data because we are running linear regressions. Some of the data might result in non-normal residuals by the nature of their data. So I applied log transformations to data that seemed related to count values or seemed log-normal (trees within 10m, density of coffee, distance to nearest tree, and wind). I applied an arcsin square-root transformation for proportion data (shade and leaves infected)

Relationships between untransformed data

chart.Correlation(roya_quadmax %>%
                    select(roya.max, coffeeDensity_m,
                           wind_m.s, shadeAvg_pc, avgTreeDist, 
                           tree_10m, evap) %>%
                    as.matrix(), histogram=TRUE, pch=19)

Relationships between transformed data

chart.Correlation(roya_quadmax %>%
                    select(roya_asqt, coffeeDens_ln,
                           wind_ln, shade_asqt, treedist_ln, 
                           tree_ln, evap) %>%
                    as.matrix(), histogram=TRUE, pch=19)

SEM models

I’ve gone ahead and fit the linear models to the transformed data since this meets the assumptions of the model. If we want in the future, we could fit generalized linear models that are appropriate for the original data distributions

SEM 1

m1 <- lm(roya_asqt ~ wind_ln + evap, data = roya_quadmax)

m2 <- lm(wind_ln ~ tree_ln, data = roya_quadmax)

m3 <- lm(evap ~ shade_asqt + wind_ln, data = roya_quadmax)

sem1 <- psem(m1, m2, m3)

summary(sem1, .progressBar = FALSE)
## 
## Structural Equation Model of sem1 
## 
## Call:
##   roya_asqt ~ wind_ln + evap
##   wind_ln ~ tree_ln
##   evap ~ shade_asqt + wind_ln
## 
##     AIC      BIC
##  45.565   76.851
## 
## ---
## Tests of directed separation:
## 
##                 Independ.Claim Test.Type  DF Crit.Value P.Value  
##      roya_asqt ~ tree_ln + ...      coef 123     1.1821  0.2395  
##           evap ~ tree_ln + ...      coef 123     2.3836  0.0187 *
##     wind_ln ~ shade_asqt + ...      coef 124    -2.0960  0.0381 *
##   roya_asqt ~ shade_asqt + ...      coef 123    -2.0270  0.0448 *
## 
## Global goodness-of-fit:
## 
##   Fisher's C = 23.565 with P-value = 0.003 and on 8 degrees of freedom
## 
## ---
## Coefficients:
## 
##    Response  Predictor Estimate Std.Error  DF Crit.Value P.Value Std.Estimate   
##   roya_asqt    wind_ln   0.0106    0.0328 124     0.3245  0.7461       0.0290   
##   roya_asqt       evap   0.0072    0.0026 124     2.7341  0.0072       0.2440 **
##     wind_ln    tree_ln   0.1336    0.1271 125     1.0508  0.2954       0.0936   
##        evap shade_asqt  -5.1995    3.3084 124    -1.5716  0.1186      -0.1378   
##        evap    wind_ln   2.5591    1.0884 124     2.3511  0.0203       0.2062  *
## 
##   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05
## 
## ---
## Individual R-squared:
## 
##    Response method R.squared
##   roya_asqt   none      0.06
##     wind_ln   none      0.01
##        evap   none      0.07
plot(sem1)

SEM 2

n1 <- lm(roya_asqt ~ evap + tree_ln, data = roya_quadmax)

n2 <- lm(evap ~ shade_asqt + tree_ln, data = roya_quadmax)

sem2 <- psem(n1, n2)

summary(sem2, .progressBar = FALSE)
## 
## Structural Equation Model of sem2 
## 
## Call:
##   roya_asqt ~ evap + tree_ln
##   evap ~ shade_asqt + tree_ln
## 
##     AIC      BIC
##  24.555   47.308
## 
## ---
## Tests of directed separation:
## 
##                 Independ.Claim Test.Type  DF Crit.Value P.Value  
##   roya_asqt ~ shade_asqt + ...      coef 123    -2.4962  0.0139 *
## 
## Global goodness-of-fit:
## 
##   Fisher's C = 8.555 with P-value = 0.014 and on 2 degrees of freedom
## 
## ---
## Coefficients:
## 
##    Response  Predictor Estimate Std.Error  DF Crit.Value P.Value Std.Estimate   
##   roya_asqt       evap   0.0068    0.0026 124     2.6293  0.0096       0.2312 **
##   roya_asqt    tree_ln   0.0554    0.0461 124     1.2024  0.2315       0.1057   
##        evap shade_asqt  -8.4275    3.3343 124    -2.5275  0.0127      -0.2234  *
##        evap    tree_ln   4.1492    1.5658 124     2.6499  0.0091       0.2342 **
## 
##   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05
## 
## ---
## Individual R-squared:
## 
##    Response method R.squared
##   roya_asqt   none      0.07
##        evap   none      0.08
plot(sem2)

SEM 3

o1 <- lm(roya_asqt ~ shade_asqt + tree_ln, data = roya_quadmax)

o2 <- lm(shade_asqt ~ tree_ln, data = roya_quadmax)

sem3 <- psem(o1, o2)

summary(sem3, .progressBar = FALSE)
## 
## Structural Equation Model of sem3 
## 
## Call:
##   roya_asqt ~ shade_asqt + tree_ln
##   shade_asqt ~ tree_ln
## 
##     AIC      BIC
##  14.000   33.909
## 
## ---
## Tests of directed separation:
## 
##  No independence claims present. Tests of directed separation not possible.
## 
## Global goodness-of-fit:
## 
##   Fisher's C = 0 with P-value = 1 and on 0 degrees of freedom
## 
## ---
## Coefficients:
## 
##     Response  Predictor Estimate Std.Error  DF Crit.Value P.Value Std.Estimate
##    roya_asqt shade_asqt  -0.2937    0.0983 124    -2.9893  0.0034      -0.2632
##    roya_asqt    tree_ln   0.1088    0.0461 124     2.3588  0.0199       0.2077
##   shade_asqt    tree_ln   0.1061    0.0409 125     2.5942  0.0106       0.2260
##     
##   **
##    *
##    *
## 
##   Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05
## 
## ---
## Individual R-squared:
## 
##     Response method R.squared
##    roya_asqt   none      0.09
##   shade_asqt   none      0.05
plot(sem3)